import java.util.ArrayList;
import java.util.Arrays;
import java.util.List;

public class PCA {
    private double[][] eigenvectors;
    private double[] eigenvalues;
    private double[] mean;

    public void fit(double[][] data) {
        int n = data.length;
        int m = data[0].length;

        // 1. 中心化数据
        mean = calculateMean(data);
        double[][] centeredData = centerData(data, mean);

        // 2. 计算协方差矩阵
        double[][] covarianceMatrix = calculateCovarianceMatrix(centeredData);

        // 3. 特征值分解
        EigenResult eigenResult = eigenDecomposition(covarianceMatrix);
        this.eigenvalues = eigenResult.eigenvalues;
        this.eigenvectors = eigenResult.eigenvectors;
    }

    public double[][] transform(double[][] data, int nComponents) {
        double[][] centeredData = centerData(data, mean);

        // 选择前nComponents个特征向量
        double[][] selectedVectors = new double[eigenvectors.length][nComponents];
        for (int i = 0; i < eigenvectors.length; i++) {
            for (int j = 0; j < nComponents; j++) {
                selectedVectors[i][j] = eigenvectors[i][j];
            }
        }

        // 投影到新的特征空间
        return matrixMultiply(centeredData, selectedVectors);
    }

    private double[] calculateMean(double[][] data) {
        int n = data.length;
        int m = data[0].length;
        double[] mean = new double[m];

        for (int j = 0; j < m; j++) {
            double sum = 0;
            for (int i = 0; i < n; i++) {
                sum += data[i][j];
            }
            mean[j] = sum / n;
        }
        return mean;
    }

    private double[][] centerData(double[][] data, double[] mean) {
        int n = data.length;
        int m = data[0].length;
        double[][] centered = new double[n][m];

        for (int i = 0; i < n; i++) {
            for (int j = 0; j < m; j++) {
                centered[i][j] = data[i][j] - mean[j];
            }
        }
        return centered;
    }

    private double[][] calculateCovarianceMatrix(double[][] data) {
        int n = data.length;
        int m = data[0].length;
        double[][] covariance = new double[m][m];

        for (int i = 0; i < m; i++) {
            for (int j = 0; j < m; j++) {
                double sum = 0;
                for (int k = 0; k < n; k++) {
                    sum += data[k][i] * data[k][j];
                }
                covariance[i][j] = sum / (n - 1);
            }
        }
        return covariance;
    }

    private EigenResult eigenDecomposition(double[][] matrix) {
        // 简化版本的特征值分解（实际应用中可使用JAMA等数学库）
        int n = matrix.length;
        double[] eigenvalues = new double[n];
        double[][] eigenvectors = new double[n][n];

        // 这里使用简化的幂迭代法（实际应用建议使用成熟的数学库）
        // 为了简化，我们假设矩阵是对称的（协方差矩阵确实对称）
        for (int i = 0; i < n; i++) {
            eigenvectors[i][i] = 1.0; // 初始化为单位矩阵
        }

        // 简化的特征值计算（实际应用应使用更稳定的算法）
        for (int i = 0; i < n; i++) {
            eigenvalues[i] = matrix[i][i];
        }

        // 对特征值和特征向量进行排序（降序）
        sortEigenvalues(eigenvalues, eigenvectors);

        return new EigenResult(eigenvalues, eigenvectors);
    }

    private void sortEigenvalues(double[] eigenvalues, double[][] eigenvectors) {
        int n = eigenvalues.length;
        for (int i = 0; i < n - 1; i++) {
            for (int j = 0; j < n - i - 1; j++) {
                if (eigenvalues[j] < eigenvalues[j + 1]) {
                    // 交换特征值
                    double tempVal = eigenvalues[j];
                    eigenvalues[j] = eigenvalues[j + 1];
                    eigenvalues[j + 1] = tempVal;

                    // 交换特征向量
                    for (int k = 0; k < n; k++) {
                        double tempVec = eigenvectors[k][j];
                        eigenvectors[k][j] = eigenvectors[k][j + 1];
                        eigenvectors[k][j + 1] = tempVec;
                    }
                }
            }
        }
    }

    private double[][] matrixMultiply(double[][] A, double[][] B) {
        int n = A.length;
        int m = B[0].length;
        int p = B.length;
        double[][] result = new double[n][m];

        for (int i = 0; i < n; i++) {
            for (int j = 0; j < m; j++) {
                double sum = 0;
                for (int k = 0; k < p; k++) {
                    sum += A[i][k] * B[k][j];
                }
                result[i][j] = sum;
            }
        }
        return result;
    }

    public double[] getExplainedVarianceRatio() {
        double total = Arrays.stream(eigenvalues).sum();
        double[] ratios = new double[eigenvalues.length];
        for (int i = 0; i < eigenvalues.length; i++) {
            ratios[i] = eigenvalues[i] / total;
        }
        return ratios;
    }

    private static class EigenResult {
        double[] eigenvalues;
        double[][] eigenvectors;

        EigenResult(double[] eigenvalues, double[][] eigenvectors) {
            this.eigenvalues = eigenvalues;
            this.eigenvectors = eigenvectors;
        }
    }
}